In this notebook, we demonstrate how to use pysofi to compresses the dynamic range of moment- and cumulant- reconstructed images using the LDRC method.
High-order cumulants or moments reconstructions result-in images with a large dynamic range of pixel intensities. This ldrc algorithm compresses the dynamic range of these reconstructions with respect to a reference image (mask) while retaining resolution enhancement.The compression is performed locally in a small window that is scanned across the image. For details of the ldrc method, please see Appendix 4 of this paper. We can pass the reference image, input image, order of reconstructions and scanning window size to the function ldrc(mask_im, input_im, order, window_size).
import sys
sys.path.append('..')
import numpy as np
import matplotlib.pyplot as plt
from functions import ldrc
from functions import reconstruction as r
from functions import visualization as v
import tifffile as tiff
%matplotlib inline
%load_ext autoreload
%autoreload 2
Here we calculate the sixth-order moments reconstruction as the input image, and the average image as the reference image. We can plot reconstructions using either matplotlib.pyplot.imshow or v.bokeh_visualization_mult.
filepath = '../sampledata'
filename = 'Block1.tif'
m6 = r.calc_moment_im(filepath, filename, order=6, frames=[0,50])
mean = r.average_image(filepath, filename)
# Plot mean and m6 with imshow
fig, axs = plt.subplots(1, 2, figsize=(10, 5))
axs[0].set_title('Reference (mean)')
axs[0].imshow(mean, cmap='pink')
axs[1].set_title('Input (6th-order moments)')
axs[1].imshow(m6, cmap='pink')
plt.show()
# Plot mean and m6 in Bokeh with hover tool and save the image as html file.
v.bokeh_visualization_mult([mean, m6],
['Reference (mean)', 'Input (6th-order moments)'],
save_option = True,
filename='Comparison of mean and m6')
As we can see from images above, the pixel values in m6 are much larger comparing to mean image, while due to the dynamic range issue, m6 appears to be dimmer and some detailos are lost.
The parameters for function ldrc.ldrc are:
mask_im: the reference image. Here we use the average image.input_im: the input image that need dynamic range compression.order: if Fourier interpolation is used when generating high-order reconstructions, this input order help to match the size between mask_im and input_im.window_size: The size of the scanning window.Please refer to the documentaion for details of each parameter for ldrc.ldrc. Window size is empirical and can be changed based on data. Here we try out different values for the window size, and choose one that yields better result.
ldrc_im = ldrc.ldrc(mask_im=mean, input_im=m6, order=6, window_size=[25, 25])
fig, axs = plt.subplots(1, 3, figsize=(15, 5))
axs[0].set_title('Reference (mean)')
axs[0].imshow(mean, cmap='pink')
axs[1].set_title('Input (6th-order moments)')
axs[1].imshow(m6, cmap='pink')
axs[2].set_title('LDRC result, window=[25,25]')
axs[2].imshow(ldrc_im, cmap='pink')
plt.show()
v.bokeh_visualization_mult([mean, m6, ldrc_im],
['Reference (mean)', 'Input (6th-order moments)', 'LDRC result, window=[25,25]'],
save_option = True,
filename='LDRC result')
ldrc_im = ldrc.ldrc(mask_im=mean, input_im=m6, order=6, window_size=[5, 5])
fig, axs = plt.subplots(1, 3, figsize=(15, 5))
axs[0].set_title('Reference (mean)')
axs[0].imshow(mean, cmap='pink')
axs[1].set_title('Input (6th-order moments)')
axs[1].imshow(m6, cmap='pink')
axs[2].set_title('LDRC result, window=[5,5]')
axs[2].imshow(ldrc_im, cmap='pink')
plt.show()
ldrc_im = ldrc.ldrc(mask_im=mean, input_im=m6, order=6, window_size=[100, 100])
fig, axs = plt.subplots(1, 3, figsize=(15, 5))
axs[0].set_title('Reference (mean)')
axs[0].imshow(mean, cmap='pink')
axs[1].set_title('Input (6th-order moments)')
axs[1].imshow(m6, cmap='pink')
axs[2].set_title('LDRC result, window=[100,100]')
axs[2].imshow(ldrc_im, cmap='pink')
plt.show()
ldrc_im = ldrc.ldrc(mask_im=mean, input_im=m6, order=6, window_size=[55, 55])
fig, axs = plt.subplots(1, 3, figsize=(15, 5))
axs[0].set_title('Reference (mean)')
axs[0].imshow(mean, cmap='pink')
axs[1].set_title('Input (6th-order moments)')
axs[1].imshow(m6, cmap='pink')
axs[2].set_title('LDRC result, window=[55,55]')
axs[2].imshow(ldrc_im, cmap='pink')
plt.show()
From images shown above, we can see that if the window size is too small, the processed image would be very niosy. If the window size is too large, the dynamic range is not compressed enough, thus the processed image is still very dim. Here we find that we get the best result when window_size=[25, 25].
You can increase the order of moments reconstruction and try out LDRC method yourself. Note that in SOFI 2.0 pipeline, if you choose to conduct Fourier interpolation as the first step, you need to change the window size accordingly based on the interpolation factor.
For cumulants-reconstructed images and moments-reconstructed images with odd orders, there appears pixels with negative values. Since the virtual fluorescence intensity is always positive, the boundaries between the negative and positive regions yield cusp-artifacts.